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o 
o 

> 

o 



O 

u 
6 

(N 
> 

00 

o 

o\ 
o 



Pratika Dayal^*, Andrea Ferrara^ & Alexandre Saro^'"* 

SISSA/International School for Advanced Studies, Via Beirut 2-4 Trieste, Italy, 34014 
^ Scuola Normale Superiore, Piazza dei Cavalieri 7, 56126 Pisa, Italy 

^ Dipartimento di Astronomia dell'Universita di Trieste, via Tiepolo 11, 1-34131 Trieste, Italy 
^ INFN, National Institute for Nuclear Physics, Trieste, Italy 



Received 2008 December 24; in original form 2008 December 24 



ABSTRACT 

We extend a previous study of Lyman Alpha Emitters (LAEs) based on hydrody- 
namical cosmological simulations, by including two physical processes important for 
LAEs: (a) Lya and continuum luminosities produced by cooling of coUisionally ex- 
cited H I in the galaxy, (b) dust formation and evolution; we follow these processes 
on a galaxy-by-galaxy basis. H i cooling on average contributes 16-18% of the Lya 
radiation produced by stars, but this value can be much higher in low mass LAEs and 
further increased if the H i is clumpy. The continuum luminosity is instead almost 
completely dominated by stellar sources. The dust content of galaxies scales with 
their stellar mass, Mdust oi M° '' and stellar metallicity, Z„, such that Mdust oc Z\-'^ . 
As a result, the massive galaxies have Lya escape fraction as low as = 0.1, with 
a LAE-averaged value decreasing with redshift: (/q) — (0.33,0.23) at z = (5.7,6.6). 
The UV continuum escape fraction shows the opposite trend with z, possibly resulting 
from dumpiness evolution. The model successfully reproduces the observed Lya and 
UV luminosity functions at different redshifts and the Lya equivalent width scatter to 
a large degree, although the observed distribution appears to be more more extended 
than the predicted one. We discuss possible reasons for such tension. 

Key words: methodsmumerical - galaxies:high redshift - luminosity function - cos- 
mology: theory 
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1 INTRODUCTION 

Lyman Alpha Emitters (LAEs) are galaxies identified by 
means of their Lya emission line, the strength and width of 
which makes the detection unambiguous to a large extent. 
Recently, there has been a massive surge in the amount of 
data available on these high redshift galaxies, specially those 
detected using the narrow band technique (Kashikawa et al. 
2006; Shimasaku et al. 2006; Dawson et al. 2007; Gronwall 
et al. 2007; Murayama et al. 2007; Ouchi et al. 2008). This 
has enabled the Lya and UV luminosity functions (LF) and 
the equivalent width (EW) distribution to be studied and 
understood statistically which has made LAEs important as 
both probes of reionization and high redshift galaxy evolu- 
tion. 

A number of theoretical models, both semi-analytic 
(Santos 2004; Dijkstra et al. 2007a,b; Kobayashi et al. 2007, 
2009; Mao et al. 2007; Dayal et al. 2008) and those involv- 
ing simulations (McQuinn et al. 2007; Nagamine et al. 2008) 
have been put forward to reproduce the observed data sets 
and put constraints on the ionization state of the IGM. How- 
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ever, the physical properties of LAEs including the star for- 
mation rates (SFR), metallicity, age and initial mass func- 
tion (IMF) have been the subject of much debate. 

To this end, in Dayal et al. (2009a) we used state 
of the art cosmological SPH simulations coupled with a 
Lya/continuum production/transmission model to infer the 
properties of LAEs. Using the SFR, age, metallicity for each 
galaxy as obtained from the simulation, a Salpeter IMF and 
the Early Reionization Model (ERM, Gallerani et al. 2007), 
we calculated the observed Lya luminosity for each galaxy 
to identify the population that would be observed as LAEs 
based on the current observational limits. We found that at 
z ~ 5.7, the SFR are in the range 2.5 — 120 M0, the ages 
range from 38 — 326 Myr, the average stellar metallicity, 
Z, ~ 0.22 Zq and the color excess, E{B -V) = 0.15. We 
thus showed that LAEs are intermediate age objects, with 
stellar metallicity higher than the generally assumed value 
of 5% Zq and dust extinction that is very much consistent 
with other works (Lai et al. 2007; Nagamine et al. 2008). 
The only free parameters left in this work were the escape 
fractions of Lya (/□,) and continuum photons (/c) from the 
galaxy, which we obtained by matching the model to the 
observations at each redshift. 
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This work extends and improves such previous study 
in many ways. The major novelty is represented by the in- 
clusion of two important features: (a) a detailed model to 
calculate dust mass evolution and the corresponding optical 
depth to continuum photons enabling us to link fc to the 
physical properties of galaxies such as the SFR, IMF, age 
and gas mass; (b) the Lya and continuum luminosities pro- 
duced by cooling of coUisionally excited H i in the galaxy 
ISM. 

The latter point is very interesting since most of the 
studies on LAEs have focused on stars as the main source of 
both Lya and continuum luminosity; an exception is Dijk- 
stra (2009), in which the case of Lya and continuum lumi- 
nosity coming from cooling of coUisionally excited H i in the 
ISM alone has been explored. Cooling of H i becomes even 
more important when one considers that Lya blobs (LABs) 
are powered either by an inflow of cooling gas (Dijkstra & 
Loeb 2009) or by supernovae (Mori, Umemura & Ferrara 
2004; Mori & Umemura 2007). 

An additional motivation of this paper is to understand 
how LAEs become progressively more dusty and how their 
dust component affects their visibility along with their ob- 
served Lya and UV continuum LPs, and line equivalent 
widths; in addition, gas emission might also have impor- 
tant effects on these quantities which need to be determined. 
The study is based on high resolution numerical simulations 
which enable us to study these effects accurately and on a 
galaxy- by- galaxy basis. 



2 SIMULATIONS 

We use state-of-the-art cosmological SPH simulations, the 
details of which can be found in Dayal et al. (2009a). In 
brief, simulations of a cosmic size {75h~^ comoving Mpc), 
have been carried out using a TreePM-SPH code GADGET- 
2 (Springel 2005) with the implementation of chemodynam- 
ics as described in Tornatore et al. (2007). The run assumes a 
metallicity-dependent radiative cooling (Sutherland & Do- 
pita 1993) and a uniform redshift-dependent Ultra Violet 
Background (UVB) produced by quasars and galaxies as 
given by Haardt & Madau (1996). The code also includes an 
effective model to describe star formation from a multi-phase 
interstellar medium (ISM) and a prescription for galactic 
winds triggered by supernova (SN) explosions, Springel & 
Hernquist (2003). In their model, star formation occurs due 
to collapse of condensed clouds embedded in an ambient 
hot gas. Stars with mass larger then 8Mq explode as su- 
pernovae and inject energy back into the ISM. These in- 
terlinked processes of star formation, cloud evaporation due 
to supernovae and cloud growth caused by cooling lead to 
self-regulated star formation. The relative number of stars 
of different mass is computed for this simulation by assum- 
ing the Salpeter (1955) IMF between 1 and IOOMq. Metals 
are produced by SNII, SNIa and intermediate and low-mass 
stars in the asymptotic giant branch (AGB). Metals and 
energy are released by stars of different masses by prop- 
erly accounting for mass-dependent lifetimes as proposed 
by Padovani & Matteucci (1993). We adopt the metallicity- 
dependent stellar yields from Woosley & Weaver (1995) and 
the yields for AGB and SNIa from van den Hoek & Groe- 
newegen (1997). Galaxies are recognized as gravitationally 
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Figure 1. Cumulative temperature distribution for interstellar 
gas averaged over all galaxies having dark matter halo masses in 
the range (in solar mass): < 10^" (solid line), IQ^^O—IOA (^sj^ort- 
dashed), iQlO-^-iO-S (dash-dotted), 10IO.8-11.2 (long-dashed), > 
10" 2 (dotted). 



bound groups of star particles by running a standard friends- 
of-friends (FOF) algorithm, decomposing each FOF group 
into a set of disjoint substructures and identifying these by 
the SUBFIND algorithm (Springel et al. 2001). After per- 
forming a gravitational unbinding procedure, only sub-halos 
with at least 20 bound particles are considered to be genuine 
structures, Saro et al. (2006). For each "bona-fide" galaxy in 
the simulation snapshot, we compute the mass- weighted age, 
the total halo/stellar/gas mass, the SFR, the mass weighted 
gas/stellar metallicity, the mass- weighted gas temperature 
and the half mass radius of the dark matter halo. 

We compute the intrinsic and observed values of the 
Lya/continuum luminosity for all the structures identified 
as galaxies in the simulation boxes at the redshifts of interest 
{z ~ 5.7, 6.6). Of these, galaxies with an observed Lya lumi- 
nosity. La ^ 10'*^'^ ergs~^ and observed equivalent width, 
EW ^ 20 A are identified as LAEs, following the current 
observational criterion for LAE identification. 

The adopted cosmological model for this work corre- 
sponds to the ACDM Universe with flm ~ 0.26, Ha ~ 
0.74, f^b = 0.0413, Us = 0.95, Ho = 73 km s^^ Mpc'^ 
and (78 = 0.8, thus consistent with the 5-year analysis of the 
WMAP data (Komatsu et al. 2009). 



3 INTRINSIC LUMINOSITIES 

As mentioned before, the huge amount of literature on LAEs 
has generally been focused on stars as the main source of 
both Lya and continuum luminosity. However, by virtue of 
using our SPH simulation, we can consistently calculate the 



^ Since we assume a simple scaling to obtain the escape fraction 
of Lya relative to continuum photons, all galaxies are classified 
as LAEs based on the EW criterion. 
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separate stellar and cooling H i gas contributions to the to- 
tal Lya and continuum luminosities, as discussed in the fol- 
lowing subsections. 



Table 1. Values of the functions Gc(T,xhi), Ga(T) and the 
corresponding intrinsic EW from cooling of coUisionally excited 
H I gas in the ISM for different halo mass ranges. 



3.1 Star powered luminosity 

Star formation produces photons with energies > 1 Ryd, 
which ionize the H i in the ISM. Due to the high density 
of the ISM, recombinations take place on a short time scale 
and this produces Lya luminosity. We represent this lumi- 
nosity as Lq, or as the Lya produced due to stars, which is 
calculated as 
2 

Ll = gQ(l - fesc)hVc, (1) 

where Q is the production rate of H i ionizing photons, fesc 
is the fraction of these ionizing photons which escape the 
galaxy without causing any ionizations, h is the Planck con- 
stant and Va is the frequency of the Lya photons. The factor 
of two-thirds arises assuming hydrogen Case B recombina- 
tion (Osterbrock 1989) and we use f^sc ~ 0.02, following the 
resuhs of Gnedin et al. (2008). 

The continuum luminosity associated with star 

formation and Q are both obtained from STARBURST99 
(Leitherer et al. 1999), using the stellar metallicity {Z), age 
(t*), IMF and the SFR (Af*) for each galaxy as obtained 
from the simulation. The use of the constant SFR is justified 
by the fact that all the galaxies identified as LAEs, at both 
the redshifts considered, have mass- weighted ages between 
30-300 Myr, which implies that recent star formation does 
not contribute significantly to the star powered Lya and 
continuum luminosity. 

As a quantitative illustration, for a galaxy with = 
200 Myr, Z = 0.2 2©, M, = IMgyr'S we find Q = 
^g53.47g-i ^j^g corresponding stellar Lya luminosity 

for f^sc = 0.02 is L; = 3.18 x lO^^ergs'^ For the same 
galaxy, the stellar continuum luminosity is L* = 3.5 x 
lO^^ergs'^A'^ which yields a stellar EW of about 91 A. 
Complete details of this calculation can be found in Dayal 
et al. (2008, 2009a). 



3.2 Luminosity from cooling H i 

In addition to the formation of stars, our numerical simula- 
tions allow us to precisely track the evolution of the amount 
of gas present in each galaxy of the studied cosmic vol- 
ume along with its temperature distribution. The cooling 
of coUisionally excited H i in the ISM of LAEs gives rise 
to both Lya and UV continuum luminosities that add to 
those powered by stars. The Lya luminosity, [erg s~^], 
produced by the cooling of neutral hydrogen gas of total 
mass Mhi, depends both on the number density of elec- 
trons and the number density of H i . We assume the hy- 
drogen gas in each galaxy to be concentrated within a ra- 
dius Tg, such that r-g = 4.5Ar2oo, where the spin parameter, 
A = 0.04 (Ferrara, Pettini & Shchekinov 2000) and r2oo is 
calculated assuming the collapsed region has an overden- 
sity of 200 times the critical density at the redshift consid- 
ered. Then, the density of electrons can be expressed as (1 — 

^ The continuum luminosity is measured in the wavelength range 
between 1250-1500A, centered at 1375 A. 
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XHi)MHi[^mH'^'J^r^\'^ and the number density of H i can 
be expressed as XHiMnillJ-'mH'i'J^rg]'^ , where fj, is the mean 
molecular weight and mn is the proton mass. Then using 
the coefficients of 7.3 x 10~"e~^^®''°°''^XH/ [ergcm^s^^] as 
given by Dijkstra (2009), the total Lya luminosity from 
cooling of coUisionally excited H i from the entire volume 
in which gas is distributed can be expressed as 

Li = 7.3 X 10" 



47rrg(/imH)^ 



where 



XHI,i 



-118400/Ti r2 
Jl,i 



(2) 



(3) 



is a nondimensional function depending on the 
H I temperature, T, and the fraction of neutral hydro- 
gen, Xhi- Further, // is the H i mass fraction mass found 
in the i-th bin of the simulated temperature distribution, 
which is discretized into A'^ = 100 bins equally spaced in 
logarithm between 10^'^ K and 10®'* K, to provide a smooth 
representation of Ga- The value of the H i fraction (xHi.i) 
in each bin is calculated by inserting the temperature (Ti) 
of that bin into the coUisional ionization-recombination 
rate equation for H i (Cen 1992). Each galaxy in the 
simulation box is first assigned to one of the five halo mass 
ranges shown in Tab. [1] Then we add its gas mass to the 
corresponding temperature bin and weigh over the total 
mass in that mass range. The values of Gc,{T,xhi) and 
galaxy-averaged temperature distributions are shown in 
Tab. [1] and Fig. [1] respectively. From the Figure, it appears 
that massive galaxies tend to have a larger fraction of their 
interstellar gas in a hot component. In eq. (2] Mhi is the 
total H I mass of the galaxy assuming a primordial gas 
composition (76% H, 24% He). 

Having determined the value of L%, the continuum lu- 
minosity, Lc [ergs~^], produced by cooling of coUisionally 
excited H i is calculated as 



LI = Li /v^)[ ^ )G.(T), 



(4) 



where we consider the 27 (2s — Is) transition to be the only 
source of the continuum and 



Gc(r) 



!:!(ls,2s) 



iV 

^ !:!(ls,2p) 

i— 1 



(5) 



Here Vc, Va are the frequencies corresponding to the contin- 
uum (1375 A) and Lya wavelengths (1216 A) respectively. 
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Figure 2. (a) The total gas mass Mgas', (b) ratio of cooling Lya luminosity from coUisionally excited H I gas, L^, to that from stars, 
L* ; (c) same ratio for the continuum luminosity at 1375 A, for z=5.7 (circles) and z=6.6 (triangles) shown as a function of the halo 
mass, Af/j, for the galaxies identified as LAEs at these redshifts. 



c is the speed of light and (jiiy /ua)du / Va is the probability 
that a photon is emitted in the frequency range v ± du/2, 
Spitzer & Greenstein (1951). Again, // is the H i mass frac- 
tion mass found in the i-th bin of the simulated temperature 
distribution, which is discretized into A'' = 100 bins equally 
spaced in logarithm between lO'^'* K and lO'' * K. We use 
<j){uc/va) = 3.1, for the continuum frequency correspond- 
ing to 1375 A. The collisional excitation rates, f2(ls, 2s) and 
Q.(ls,2p) are calculated for each temperature bin by using 
the results of Anderson et al. (2000). Again using the tem- 
perature distribution for all galaxies as shown in Fig. [T] the 
mass-weighted value of this ratio varies with the halo mass 
range as shown in Tab. [1] However, as will be shown later, 
the continuum from the stars always dominates over L®, 
which makes our results insensitive to the exact values of the 
collisional excitation rates and continuum wavelength cho- 
sen 0. Using these conversion factors, the calculated value 
of intrinsic EW of the Lya line from gas alone has a value 
between 887-1304 A, depending on the halo mass, as shown 
in Tab. [T] 

Once that La and Lc from both the stars and ISM gas 
component have been determined, the total intrinsic equiv- 
alent width, EW^"^, can be calculated as 

EW'"' = 4^-^. (6) 

The observed equivalent width, EW , is affected by dust in 
the ISM (modifying both the continuum and Lya luminos- 
ity) and transmission through the IGM, which only damps 
the LyQ luminosity. The observed EW in the rest frame of 
the galaxy is calculated as 

EW = EW""' [ ^]Ta, (7) 



The intrinsic EW from cooling of coUisionally gas alone de- 
creases from 1449 to 681 A, as the continuum wavelength in- 
creases from 1250 to 1500 A. However, since the stellar contin- 
uum always dominates that from the gas, we show the results at 
a wavelength corresponding to 1375 A, for convenience. 



where fa/fc represents the differential effect of dust on 
Lya and continuum photons; a fraction Tq of the Lya lu- 
minosity emerging from the galaxy is transmitted through 
the IGM. The transmission through the IGM depends on a 
number of factors including the SFR, the escape fraction of 
H I ionizing photons, the size of the Stromgren sphere built 
by the galaxy, the global value of the fraction of H i and the 
extent to which sources are clustered. Complete details of 
the calculation of Ta can be found in Dayal et al. (2009a). 
The effect of dust on both the continuum and Lya luminos- 
ity is discussed in Sec. |4l 

We find that the contribution of to the total intrin- 
sic Lya luminosity is not negligible; the average value of 
L^/La — (0.16,0.18) at 2 = (5.7,6.6) respectively, as seen 
from Panel (b). Fig. [2] The ratio L^/L'^, on the other hand 
is about a factor ten smaller, i.e. (0.01,0.02) at z = (5.7, 6.6), 
respectively (Panel (c) of same Figure). Adding and L^, 
therefore, has the effect of boosting up the intrinsic EW. As 
shown in Dayal et al. (2009a) , the average intrinsic EW pro- 
duced by stars only at z ~ 5.7 — 95 A. However, as shown in 
Tab. [2] we find that adding boosts the average intrinsic 
EW to 106 A, i.e., by a factor of about 1.12 as expected. 

It is interesting to see that even though the gas mass 
(Panel (a) , Fig. [2| is almost similar for a given halo mass 
at 2 ~ 5.7 and 6.6, both and are slightly larger at 
z ~ 6.6. This can be explained noting that the virial ra- 
dius decreases with increasing redshift, which leads to a cor- 
responding increase in the H i density. The higher density 
leads to an increase in the collisional excitation rate and 
hence, the cooling rate, producing more and L^. 

Further, as shown in Dayal et al. (2009a), the SFR rises 
steeply with increasing halo masses, which leads to L* and 
Lg becoming increasingly dominant as compared to La,Lc. 
Hence, the contribution of ISM gas towards both the Lya 
and continuum luminosity decreases with increasing halo 
mass. The dip in (and the corresponding dip in Lf ) for 
halo masses between lO^"'^ and 10^^°'* Mq is caused by the 
lower value of Gc,{T,xhi) for these masses, as seen from 
Tab.[Tl 
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4 DUST MODEL 

Dust is produced by supernovae and evolved stars in a 
galaxy. However, several authors (Todini & Ferrara 2001; 
Dwek et al. 2007) have shown that the contribution of AGB 
stars becomes progressively less important and at some point 
negligible towards higher redshifts {z ^ 5.7). This is because 
the typical evolutionary time-scale of these stars (S5 1 Gyr) 
becomes longer than the age of the Universe above that red- 
shift. However, under certain conditions thought to hold in 
quasars, in which extremely massive starbursts occur, the 
contribution of AGB can become important somewhat ear- 
lier (see Valiante et al. 2009). We therefore make the hy- 
pothesis that the dust seen in LAEs at z J5 5.7 is produced 
solely by type II supernovae. We compute the amount of 
dust in each LAE by applying the evolutionary model de- 
scribed below by post-processing the simulation outputs. 

As mentioned before, for each simulated galaxy the 
SFR, M* , the total stellar mass, M* and the mass- weighted 
stellar age, t, are available. We then use the average SFR 
over the SF history of the galaxy as (M*) = M^/tt to cal- 
culate the dust mass. This is because the SNII rate, and 
hence the dust content of the galaxy depend on the entire 
SFR history and not just the final values obtained from the 
simulation. The SNII rate is estimated to be 7(M*), where 
7 ~ (54M0)~^ for a Salpeter IMF between the lower and 
upper mass limits of 1 and 100 respectively (Ferrara, 
Pettini and Shchekinov 2000). We have assumed that the 
progenitors of SNII have a mass larger than 8 and in or- 
der to simplify the calculation, we have adopted the instan- 
taneous recycling approximation (IRA), i.e., the lifetime of 
stars with mass larger than 8 M0 is zero. 

Using this SNII rate, we can calculate the evolution of 
the total dust mass, Mdustif), in the galaxy as 

dM^ . ,.7(M.) - ^ - ^(M,), (8) 
at rdest(t) Mg(t) 

where the first, second and third terms on the RES repre- 
sent the rates of dust production, destruction and astra- 
tion (assimilation of a homogeneous mixture of gas and 
dust into stars) respectively; the initial dust mass is as- 
sumed to be zero. Further, yd is the dust yield per SNII, 
Tdest is the timescale of dust destruction due to SNII blast 
waves and Mg{t) is the gas mass in the galaxy at time t. 
For yd, the average dust mass produced per SNII, we adopt 
a value of 0.5 M© (Todini & Ferrara 2001; Nozawa et al. 
2003, 2007; Bianchi & Schneider 2007). The amount of gas 
left in the galaxy, Mg{t), at any time t is calculated as 
dMg{t)/dt — — (M,) and each galax;y is assumed to have 
an initial gas mass given by the cosmological baryon to dark 
matter ratio Qb/ilm- The final gas mass so obtained for each 
galaxy is very consistent with the value obtained from the 
simulation snapshot. This is primarily because LAEs are 
much larger than dwarfs and hence, are not expected to 
have very strong outfiows. This also justifies neglecting the 
outflow term due to galactic winds in eq. ([8]). 

The destruction timescale, Tdest{t) is estimated follow- 
ing the resuhs of McKee (1989), Lisenfeld & Ferrara (1998) 
as: 

Mgit) 



Tdest (t) = 



(9) 



7 (M.>eMs( 100 kms-i) 
Here, e is the efficiency of dust destruction in a SN-shocked 
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Figure 3. (a) Dust mass and (b) fa as a function of halo mass, 
Mfi, for z=5.7 (circles) and z=6.6 (triangles) for the galaxies iden- 
tified as LAEs at these redshifts using p = (1.5, 0.6) for all galaxies 
at 2 ^ (5.7,6.6). For details, refer to text in Sec. |4] 



ISM, for which we adopt value ~ 0.4; this is a reasonable 
estimate between the values of 0.1 and 0.5 found by Mc- 
Kee (1989) and Seab & ShuU (1983) for varying densities of 
the ISM and magnetic field strengths. Ms(100kms~^) is the 
mass accelerated to 100 km s~^ by the SN blast wave and 
has a value of 6.8 x 10^ M© (Lisenfeld & Ferrara 1998). 

Since we consider a scenario in which the dust amount is 
regulated only by the SFR, which scales with the halo mass, 
the dust amounts too scale in the same way, as seen from 
Panel (a) of Fig. |3] However, while the dust mass increases 
linearly with the halo mass for most of the objects, the scat- 
ter at the low mass end arises from feedback regulation of 
the star formation in these objects. 

Once the final dust mass, Mdust{t,), is calculated for 
each galaxy in the simulation, we can translate this into an 
optical depth, Tc, for continuum photons as 

_ 3Ed 



(10) 



where Ed is the dust surface mass density, a and s are the ra- 
dius and material density of graphite/carbonaceous grains, 
respectively {a — 0.05fim, s = 2.25 gcm""^; Todini & Fer- 
rara 2001; Nozawa 2003). The dust surface mass density is 
calculated assuming that grains are concentrated in a ra- 
dius proportional to the stellar distribution, r^, following 
the results of Bolton et al. (2008) , who have derived fitting 
formulae to their observations of massive, early type galaxies 
between 2 = 0.06 — 0. Though it is not an entirely robust 

The relation used to obtain the stellar distribution radius, re, 
from the observed V band luminosity, Ly is 



log re 

{Mdiru/lO^'-MQ) 



1.28 log (Te2 - 0.77 log /e 

0.691(Lv/10"LQ)l■2^ 



0.09, 



(11) 
(12) 
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estimate, we extend this result to galaxies at z 5.7, 6.6 due 
to the lack of observations about the stellar distribution in 
high-redshift galaxies. 

We calculate — Mdust{t*) / {t^t"^) , where rd is the 
effective radius of dust. The best fit to the observed UV 
LFs requires the dust distribution scale to be quite similar to 
that of the stellar distribution, such that ra = (0.6, 1.0) x re 
at z - (5.7, 6.6). 

This optical depth can be easily converted into a value 
for the escape fraction of continuum photons assuming a 
slab-like distribution of dust, such that 



We are hence able to determine fc for each individual galaxy, 
depending on its intrinsic properties, i.e. the SFR, metallic- 
ity, age, gas mass and IMF. 

As /c is fully determined from the above procedure, we 
are left with a single free model parameter, fa, to match 
the Lya LF obtained from our calculations to the data. The 
relation between fa and fc can be written as 

fa=p{A,,,C)fc, (15) 

where p depends not only on the adopted extinction curve, 
A\ , but also on the differential radiative transfer effects act- 
ing on both Lyof and UV continuum photons in an inho- 
mogeneous medium described by a clumping factor C = 
(ri^)/{n)^. While for the former dependence, some pieces of 
evidence exist that the SN dust we assume here can suc- 
cessfully be used to interpret the observed properties of the 
most distant quasars (Maiolino et al. 2006) and gamma- 
ray bursts (Stratta et al. 2007), the way in which inhomo- 
geneities differentially enhance Lya with respect to UV con- 
tinuum luminosities through the so-called "Neufeld effect" 
(Neufeld 1991; Hansen & Oh 2006; Verhamme et al. 2008; 
Finkelstein et al. 2009; Kobayashi et al. 2009) is still highly 
debated. A report of a decreasing trend of fa with E{B — V) 
has been recently published (Atek et al. 2009); our results 
nicely match that empirical relation, as shown in Fig. [l] 
We therefore tentatively assume C = 1 and compute p sim- 
ply from our SN extinction curve, which corresponds to a 
value of Rv = 2.06, obtaining p{A\, C = 1) — 0.8, indepen- 
dently of the LAE mass or luminosity. Although the best fit 
to the observed Lya LF favors values that vary around 0.8 
(p = 1.5,0.6 at z = 5.7,6.613, it is unclear if these devia- 
tions can be interpreted as a genuine imprint of an evolving 
clumpy ISM. 

We now return to the results obtained using the dust 

<Te2 = (2GM<i,„/re)°-""Ml00km/s], (13) 
/e = Lv/(27rr2)[lo9Lv/(i0kpc2)]. (14) 

Here, Ly (the luminosity between 4644 and 6296 A in the ob- 
server's frame) is obtained using STARBURST99 for the SFR, 
age, metallicity and IMF of the galaxy under consideration. For 
a galaxy at z ~ 5.7 which has t, = 200 Myr, Z = 0.2 Z©, 
M* = lM0yr~i, the V band luminosity is 5.05X lO^^ergs"^ 
and the corresponding re = 0.35 kpc. Although we have used the 
V band luminosity in the observer's frame in this work, using the 
rest frame V band luminosity does not affect our results in any 
way, using the same values of the free parameters. 
^ We set p = 1 when fa ^ 




0^2 0.4 0.6 0.8 



E(B-V) 

Figure 4. The escape fraction of Lya photons, fa as a function 
of the dust color excess, E{B — V). Filled points (circles, trian- 
gles, squares) show the observed values for Lya galaxies from the 
GALEX LAE sample (Deharveng et al. 2008), lUE local starburst 
sample (McQuade et al. 1995; Storchi-Bergmann et al. 1995) and 
Atek et al. (2008) sample respectively. The dashed line is the best 
fit to all the observations (see Atek et al. 2009), the dot dashed 
line shows the best fit from Verhamme et al. (2008) obtained from 
the spectral fitting of 2 ~ 3 LBGs and the dotted line is obtained 
following the Seaton law (1979). The relations obtained from this 
work using p = (1.5, 0.6) for LAEs are shown with empty circles 
and triangles for z ~ (5.7, 6.6), respectively. Details in Sec. |4] 




4 5 6 7 

Log {M.^JM^} 

Figure 5. Relation between (a) stellar mass, M*, (b) stellar 
metallicity, Z» and dust mass, M^ustt for z=5.7 (circles) and 
z=6.6 (triangles) for the galaxies identified as LAEs at these red- 
shifts. 
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model. For galaxies identified as LAEs at z ~ 5.7, the es- 
cape fraction of Lya photons, fa, decreases by a factor of 
10, from 1 to about 0.1 as the halo mass increases from 10^ ° 
to lO'^^'^ M0 as shown in Panel (b) of Fig.|3] consistent with 
the results obtained by Laursen et al. (2009). However, for 
the few galaxies between 10^^ '^ and 10^^ * Mq, /q increases 
again. This is because the optical depth, r oc Mdust/r1, 
where scales with SFR (see eq. 11). Since the SFR rises 
steeply with halo mass for Mj, ^ lO^'^'^ Mq, the optical 
depth decreases, even though the total dust mass increases. 
The decreased optical depth then leads to the larger escape 
fraction of continuum and hence, Lya photons seen. 

We find that the galaxies with the largest stellar mass 
are more dust rich, as shown in Panel (a) of Fig. [5] for 
the bulk of the galaxies, the dust mass, Mduat is related 
to the stellar mass, M, , as Mdust tx M°'^. However, the low 
stellar mass end shows a large scatter in the dust content, 
likely arising from feedback regulation leading to a scatter 
in the star formation rates for smaller halos/galaxies. As 
expected, galaxies with the highest stellar metallicities have 
the largest dust masses; the dust mass is related to the stellar 
metallicity, Z,, as Mduat oc as shown in Panel (b) of 
Fig. [5l Again, the scatter for the low metallicity galaxies is 
explained by the corresponding scatter in the star formation 
rates at this end. 

We find that {fc) ~ (0.23, 0.38) at 2 ~ 5.7, 6.6 respec- 
tively. This is consistent with the values of (0.22, 0.37) ob- 
tained in Dayal et al. (2009a). The Lya contribution from 
the gas allows smaller galaxies to become visible as LAEs: 
the minimum halo mass of LAEs shifts from (10^°-^)Mo 



to (10« ^ 10"-^) at 



(5.7,6.6). Since smaller galaxies 



{Mh ^ 10 " M0) have lower SFR and dust content, they 
tend to have a higher fc value. However, only about 10% 
of the galaxies we identify as LAEs in the simulation have 
Mh < lO^^'^M© and hence, they are too few to vary the 
mean /c. 



5 OBSERVATIONAL IMPLICATIONS 

As mentioned before, cooling of collisionally excited H i can 
produce a significant amount of Lya luminosity in the ISM. 
However, since galaxies become more compact and denser 
towards high redshift, this leads to a corresponding increase 
in L%. This can be clearly seen from the cumulative Lya 
LF (Fig. [S} where ignoring the gas contribution leads to a 
slightly larger under estimation of the LF at 2 ~ 6.6 com- 
pared to 2 ~ 5.7. 

As mentioned in Sec. [S] cooling of collisionally excited 
H I does not contribute much to the continuum luminosity, 
with the average contribution to the total continuum lumi- 
nosity being ^ 10 % at either of the redshifts considered. 
Adding Lf, therefore, produces a negligible effect on the 
UV LF, as seen clearly from Fig.[71 where the UV LFs, both 
including and excluding the gas contribution overlap. 

The free parameters in this work are the the radius of 
the dust distribution, rd relative to the stellar distribution 
radius, re and the escape fraction of Lya photons relative 
to continuum photons. Once that these two parameters are 
obtained by calibrating our results to the UV and Lya LFs, 
we have no more free parameters left. Using the values of 
fa/ fc, Ta and the intrinsic EW for each galaxy we identify 
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z~6.6 
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42.5 43 

Log {Ljerg s-') 

Figure 6. Cumulative Lya LF for the ERM. The panels are 
for z ~ 5.7 and 6.6. Points represent the data at two differ- 
ent redshifts: z ~ 5.7 (Shimasaku et al. 2006) (circles), z ^ 6.6 
(Kashikawa et al. 2006) with downward (upward) triangles show- 
ing the upper (lower) limits. Solid (dashed) Lines refer to model 
predictions at 2: ~ 5.7, 6.6 for the parameter values in text in- 
cluding (excluding) the contribution from cooling of collisionally 
excited H I in the ISM. Shaded regions in both panels show pois- 
sonian errors. 
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Figure 7. UV LAE LF for the ERM. Points represent the data at 
two different redshifts: z ~ 5.7 (Shimasaku et al. 2006) (circles), 
z ~ 6.6 (Kashikawa et al. 2006) (triangles). Solid (dashed) lines 
refer to model predictions including (excluding) the contribution 
from cooling H I in the ISM at the redshifts (from top to bottom): 
z ~ 5.7, 6.6, for the parameter values in text. The vertical dashed 
(dotted) lines represent the observational 2a (Scr) limiting mag- 
nitude for z = 5.7 {z = 6.6). The shaded regions in both panels 
show the poissonian errors. 
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Table 2. For all the LAEs comprising the Lya LF at the redshifts shown (col. 1), we show the range of halo mass (col. 2), the range 
of mass weighted ages (col. 3), the range of SFR (col. 4), the range of gas mass (col. 5), the average color excess (col. 6), the average 
transmission of Lya photons through the IGM (col. 7), the average intrinsic EW (col. 8), the average EW of the line emerging from the 
galaxy (col. 9) and the average value of the observed EW (col. 10). 



z 


Mh 

[Mq] 


[Myr] 


[MQj/r-i] 




{E{B-V)) 




[A] 


[A] 


{EW) 

[A] 


5.7 
6.6 


]^q9.0-11.8 
]^q9.4-11.6 


44.0 - 326.2 
23.3 - 218.1 


0.8 - 120 
1.6 - 46.4 


]^q9. 0-11.0 
]^q9.1-10.7 


0.14 
0.09 


0.48 
0.48 


106.0 
117.2 


157.2 
70.3 


75.1 
34.6 



3 r^-^ 




42.5 43 43.5 44 

Log (L„/erg s-') 



Figure 8. Observed EWs from Shimasaku et al. (2006) (circles) 
and model values of the observed EWs (astrexes) as a function of 
the observed Lya luminosity. 

as a LAE, we obtain the emergent and observed EWs at 
z ~ 5.7. Due to the cooling of collisionally excited H i , 
the total Lya luminosity increases by a factor of about 1.12 
while the continuum almost remains unchanged. This leads 
to intrinsic EWs in the range 88.8-188.8 A for the galaxies 
we identify as LAEs. The average value of the intrinsic EW 
at 2 ~ 5.7 also increases to 106.0 A compared to 95 A that 
we found in Dayal et al. (2009a), not including gas emission. 

The EWs emerging from the galaxy at z ~ 5.7 are fur- 
ther enhanced and have a mean value of 157.2 A because 
of the larger escape fraction of Lya photons with respect to 
continuum photons {fa/fc = 1.5). However, since the Lya 
luminosity is attenuated by transmission through the IGM, 
with {Ta) ~ 0.48, the observed EWs are finally reduced to 
the range 44.4-106.9 A; with an average value of 75.1 A. 

We show the EWs as a function of the observed Lya 
luminosity in Fig. [S] Though the numerical data is concen- 
trated at Log(i5VK/A) ~ 1.9, it shows some scatter, mostly 
induced by the different physical properties of galaxies as age 
and metallicity. However, we are unable to reproduce some 
of the largest and smallest EW values. This could be due to 
some of the assumptions made in the modelling of in- 
cluding: (a) we have calculated for a homogeneous ISM. 
This assumption seems to represent a fairly good description 



of the \ow-z data discussed in Fig. 4. However, we cannot 
exclude the possible presence of gas inhomogeneities which 
would increase the coUisional excitation rates and hence in- 
crease Lq. Since we are investigating a large cosmological 
volume, resolving the small scale structure of the ISM is be- 
yond the possibilities of the present work, and our results in 
this sense should be considered as a lower limit to the true 
emissivity of the gas. (b) depends on the radius of gas 
distribution and the simple relation we have assumed might 
be different for these high redshift galaxies and (c) we have 
assumed a very simple scaling between and fc. In real- 
ity, this should depend on the geometry and distribution of 
dust within the ISM and this would induce scatter in the 
EW distribution. 



6 SUMMARY 

Coupling SPH simulations with a Lya produc- 
tion/transmission model, in Dayal et al. (2009a) we 
studied the evolution in the properties of Lya emitters 
at high redshifts. In this paper, we further extend that 
work by including two important ingredients, based on the 
intrinsic properties of each galaxy: (a) Lya and continuum 
luminosities produced by cooling of collisionally excited 
H I with electrons in the ISM and (b) a dust model to 
calculate the optical depth and hence, the escape fraction 
of continuum photons (and color excess) from each galaxy. 

We find that the Lya luminosity from cooling ISM gas 
(Lq) constitutes a non- negligible addition to the value from 
stars (I/a); the average value of L^/L^ ~ (0.16,0.18) at 
2 ~ (5.7, 6.6). However, the continuum from stellar luminos- 
ity (L*) always dominates over that from cooling gas (L^); 
the average value of Lg/L* ~ (0.01,0.02) at 2 - (5.7,6.6). 
The contribution of the cooling of H i to both the Lya and 
continuum luminosities increases with redshift because of 
the increasing H i gas density, which leads to a larger exci- 
tation and hence, cooling rate. Since the SFR rises steeply 
with increasing halo mass, L* and L* also increase faster 
as compared to and and hence dominate in massive 
galaxies. 

The dust budget of a given LAE is determined by sev- 
eral physical processes (dust formation, destruction and as- 
tration) mainly regulated by the evolution of massive stars. 
Hence the dust mass increases with increasing halo mass, 
following the behavior of the star formation rate. Galaxies 
with the largest stellar mass, M* have also the largest dust 
content, with a scaling relation given by Mdust oc The 
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same trend of increasing dust content is also observed in 
terms of the stellar metallicity, Mdust oc Zl''' . 

We calculate the escape fraction of continuum pho- 
tons using a slab-like dust distribution model and 
graphite/carbonaceous grains (radius a = 0.05/im, density 
p = 2.25gcm~^). We then assume a mass-independent re- 
lation to go from /c to fa such that fa = (1.5, 0. 6) /c at 
z ~ (5.7, 6.6) respectively. The increase in the dust mass 
with increasing halo mass leads to a decrease in the escape 
fraction of Lya photons [fa) from the galaxy, with the value 
dropping from ~ 1 to 0.1 as the halo mass increases from 
10^-° to 10"-^ M0 at 2 ~ 5.7. 

Using a SN dust extinction curve (Todini & Ferrara 
2001), fa = 0.8/c. This means that the ISM dust is quite 
homogeneous for fa ^ 0.8/c; a situation we find at z ~ 6.6. 
However, no extinction curve can produce the value of fa > 
fc, which hints at the fact that the dust could be clumped 
at z ~ 5.7 as proposed by Neufeld (1991) and Hansen & Oh 
(2006). 

Using two free parameters to translate fc to fa and re to 
rd, we are able to reproduce the UV and Lya LPs reasonably 
at both the redshifts considered. With no more free parame- 
ters, we calculate the intrinsic, emergent and observed EWs 
from the galaxies we identify as LAEs. Adding Lya from the 
cooling H I , the average intrinsic EW at z ~ 5.7 increases 
by a factor of about 1.12 and has a mean vale of 106.0 A. 
The emergent EW is further enhanced to an average value of 
about 157.2 A at 2 ~ 5.7 since fa/fc ~ 1.5 at this redshift. 
However, attenuation of the Lya luminosity H i in the IGM 
reduces the mean observed EW to 75.1 A. 

The observed EWs (Shimasaku et al. 2006) show a large 
scatter when plotted against the observed Lya luminos- 
ity. By virtue of using a Lya production/transmission and 
dust model that depends on the intrinsic properties of each 
galaxy, we also obtain a scatter in the same plot. 

We conclude by discussing some unsolved points still 
left in our calculations. Though we are able to reproduce 
the Lya LP at 2: ~ 5.7 to a reasonable extent, we are still 
unable to capture the shape of the slope completely. Purther, 
though we get a scatter in the EWs, we are still unable to 
produce the largest and smallest EWs. These could be due 
to the simplifying assumptions adopted in our modelling. 
These include: (a) a simple scaling relation between the gas 
distribution scale and the virial radius for each galaxy. How- 
ever, a spread in this relation would vary the calculated; 

(b) assuming a clumpy ISM would increase La as compared 
to the value we get assuming a homogeneous distribution; 

(c) a simple scaling between fa and fc- Note, however, that 
a relation dependent on the individual ISM conditions of 
each galaxy would induce a scatter in the EW distribution; 

(d) inflows/outflows, could vary the EWs and probably also 
leave an imprint on the Lya LF and (d) a detailed radia- 
tive transfer model of reionization which could affect our 
calculation of Ta ■ All these issues will be investigated in fu- 
ture works with the aim of building a complete and coherent 
model for LAEs. 
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